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DOMAIN  DECOMPOSITION  ALGORITHMS  FOR  THE  BIHARMONIC 

DIRICHLET  PROBLEM 

XUEJUN  ZHANG' 

Abstract.  We  consider  additive  Schwarz  methods  for  the  biharmonic  Dirichlet  problem  and  show  that 
the  algorithms  have  optimal  convergence  properties  for  some  conforming  finite  elements.  Some  multilevel 
methods  are  also  discussed. 

1 .   Introduction.  We  are  interested  in  solving  the  following  biharmonic  Dirichlet  prob- 
lem in  a  plane  region 


(1) 


f   A^u     =     /  in  n, 

u    =    go        on        dQ, 

^  =  91      on      an. 


It  is  convenient  to  work  with  the  weak  formulation:  Find  u  G  Hq{Q.)  such  that 

(2)  a{u,v)^f{v),     \/veH^{il), 

where  /  is  a  bounded  linear  functional  on  H^(^)  and  a{u,v)  is  a  symmetric,  continuous, 
/To-elliptic  bilinear  form.  Two  examples  of  such  bilinear  forms  are 

(3)  a{u,v)=   /  AuAv  dx, 
and 

(4)  aiu,v)  =  JjAuAv  +  (1  -  .)(2^^^^^^  -  ^^  -  ^^)}  dx, 

where  0  <  a  <  1/2  is  Poisson's  coefficient  of  the  plate.  The  first  one  arises  in  Fluid 
Dynamics,  and  the  second  provides  a  variational  formulation  of  the  Clamped  Plate  Problem. 
The  rest  of  the  paper  is  organized  as  follows.  In  section  2,  we  introduce  some  standard 
conforming  finite  element  approximations.  In  section  3,  we  study  some  additive  Schwarz 
algorithms  and  establish  optimal  convergence  properties  of  the  algorithms.  In  section  4,  we 
study  a  multilevel  algorithm  for  the  biharmonic  equation.  In  this  brief  paper,  we  do  not 
provide  proofs  and  details  of  the  algorithms;  cf.  Zhang  [11]  for  further  discussion  of  the 
algorithms  and  some  numerical  experiments. 

2.   Finite  Elements  for  the  Biharmonic  Equation. 

2.1.  Formulation.  We  triangulate  the  domain  Q.  into  non-overlapping  regions  called 
elements,  generally  triangles  or  rectangles.  Let  V^  be  a  space  of  piecewise  polynomials  with 
respect  to  the  triangulation.  The  finite  element  solution  Uh  G  V^  satisfies 

(5)  a{uH,(l>h)  =  f{4>h),     V</.;,GV^ 


•    Courant  Institute,  251  Mercer  St,  New  York,  NY  10012.    xuejun@widlund.nyu.edu.    This  work  was 
supported  in  part  by  the  National  Science  Foundation  under  Grant  NSF-CCR-8903003. 


Let  {(f)i}  be  the  nodal  basis  for  V''.  Then  Uh  can  be  represented  as 

i 

Thus,  we  obtain  a  linear  system  for  x,  the  degrees  of  freedom  of  u^, 

K^x  =  b, 

where  A'''  =  {a( </>,,(/>_,)}  and  6,  =  /(</>,). 

The  stiffness  matrix  K^  is  symmetric,  positive  definite.  After  a  proper  scaling,  its 
condition  number  k{K'^)  =  0{h~'*).  Since  the  system  are  usually  very  large,  and  the 
condition  number  of  A"''  is  also  very  large,  solving  the  system  can  be  very  expensive.  Many 
preconditioners  have  been  designed  for  A'''.  Among  them,  the  additive  Schwarz  methods 
studied  in  this  paper  seem  to  be  particularly  successful  and  promising. 

2.2.  Some  Conforming  Elements.  For  biharmonic  equation,  the  finite  elements  are 
all  relatively  complicated.  In  this  paper,  we  restrict  ourselves  to  some  standard  conforming 
elements.  In  particular,  we  consider  the  Argyris  triangle  V^,  the  Bell  triangle  Vg  and  the 
bicubic  element  VA.  These  elements  are  complicated  but  among  the  simplest  conforming 
elements  for  the  biharmonic  equation;  cf.  Ciarlet[5]. 

The  Argyris  element  consists  of  continuous  differentiable  functions,  the  restriction  of 
which  to  any  element  is  in  V5.  The  degrees  of  freedom  for  the  Argyris  element  in  a  triangle 
with  vertices  a,,  i  =  1,2,3,  are  given  by 

(     £:p(a.)M<2,    £-M6.)     ), 

where  6,  is  the  midpoint  of  the  edge  0^0^,  and  ra,  is  the  outward  normal  of  "ojOk.    The 
number  of  the  degrees  of  freedom  for  one  triangle  is  21. 

It  is  easy  to  see  that,  in  general,  the  normal  derivatives  of  an  Argyris  element  is  a 
polynomial  of  degree  4.  Let  Vb  denote  the  subspace  of  V5  formed  by  those  polynomicds  of 
V5  whose  normal  derivatives  along  each  side  of  a  triangle  are  polynomials  of  degree  3int,t 
being  the  abscissa  along  an  axis  containing  the  side.  We  note  that  P4  C  Vb  C  Vs-  The  Bell 
element  consists  of  C^  functions  whose  restrictions  to  a  triangle  are  in  Vb-  The  degrees  of 
freedom  for  the  Bell  element  are  given  by 

{  ^p(a.),H<2}. 

K  the  domain  is  built  from  rectangles,  then  we  can  also  use  the  bicubic  element,  known  as 
the  Bogner-Fox-Schmit  rectangle  in  the  engineering  literature.  It  is  a  space  of  C*  functions 
whose  restrictions  to  a  rectangular  element  are  in  Q3  =  span{x'y^, 0  <  i,j  <  3}.  The 
degrees  of  freedom  of  bicubic  element  are  given  by 

,     t    ^     dp  dp  _^P_,     X  -, 

{  Pia,),  ^(a.),  ^(a.),  ^^;^(a.)  }• 

3.  Additive  Schwarz  Methods  for  the  Biharmonic  Problem.  In  this  section,  we 
study  the  additive  Schwarz  method  for  the  biharmonic  problem.  We  consider  the  bicubic 
element,  the  Argyris  element  and  the  Bell  element. 


The  additive  Schwarz  schemes  were  designed  by  Dryja  and  Widlund  [6]  and  Matsokin 
and  Nepomnyaschikh  [9].  The  optimal  convergence  properties  of  the  algorithm  were  es- 
tablished for  second  order  self-adjoint  elliptic  problems,  see  Dryja  and  Widlund  [6,7,8]. 
Generalizations  to  the  nonsymmetric  or  indefinite  cases  have  been  made  by  Cai  and  Wid- 
lund; cf.  [2,3,4].  We  show  in  this  paper  that  the  condition  number  of  certain  additive 
Schwarz  methods  for  the  biharmonic  equation  is  uniformly  bounded. 

Suppose  that  the  finite  element  space  V  can  be  written  as  a  sum  of  subspaces. 

Instead  of  solving  the  original  finite  element  equation,  in  the  additive  Schwarz  scheme,  we 
solve 

Puh  =  (Pvo  +  Pv,  +  •  •  •  +  Pvr,)uh  =  9k 

for  some  5/,.  Here  Py,  :  V  — >  Vi,  is  a  projection  defined  by 

(6)  aiPv,u,<t>)  =  a{u,<f>),     V<^  G  V.. 

The  natural  question  is  how  to  find  decompositions  of  V'',  and  what  properties  of  the 
decomposition  give  optimal  algorithms. 

As  for  the  second  order  cases,  a  coarse  problem  is  crucial  in  the  algorithms.  In  the 
second  order  case,  an  obvious  candidate  for  the  coarse  subspace  is  the  space  V"  associated 
with  a  coarse  triangulation  T^ .  However,  for  the  biharmonic  case,  when  the  Argyris  and 
Bell  elements  are  used,  the  coarse  finite  element  space  V^  (V^)  is  not  a  subspace  of  V^ 
(V^).  Therefore,  we  cannot  use  V^  or  V^  as  our  coarse  subspaces  and  a  new  coarse 
subspace  has  to  be  found. 

In  the  case  of  the  iterative  substructuring  methods,  the  situation  is  even  worse.  We 
recall  that  these  are  domain  decomposition  algorithms  which  use  nonoverlapping  subregions; 
cf.  Widlund  [7].  When  we  use  the  bicubic  element  and  thus  V^  C  V^  can  be  used  as 
coarse  subspace,  the  direct  generalization  of  certain  algorithms  designed  for  second  order 
problems  results  in  algorithms  with  condition  numbers  which  grow  at  least  like  (H/h)'^. 
Better  algorithms  are  obtained  by  adding  certain  vertex  spaces  to  the  space-decomposition. 
We  will  not  discuss  the  iterative  substructuring  method  further  in  this  paper. 

The  difficulty  of  proving  the  optimality  of  the  algorithms  are  due  to  the  presence  of  the 
high  order  derivatives  in  the  elements.  The  tools  that  work  for  second  order  equations  and 
linear  element  cannot  be  used  here. 

3.1.  The  bicubic  Element.  We  first  triangulate  the  domain  il  into  nonoverlapping 
rectangles  fi,,  i  =  1,---,N,  to  obtain  a  coarse  triangulation  T^  =  {fit}^.  Then  each 
rectangle  Q,  is  further  divided  into  smaller  rectangles  r,  to  obtain  the  fine  triangulation 

r'  =  {r,}. 

We  assume  that  fi  can  be  decomposed  into  overlapping  subdomains  il  =  U.-jfi,  and 
the  decomposition  satisfies 

Assumption  3.1.   The  decomposition  Q.  =  U^jfii  satisfies 

•  ddi  aligns  with  the  boundaries  of  fine  elements,  i.e.  ft,  is  the  union  of  some  ele- 
ments Tj 

•  {ft}i^i  forms  a  finite  covering  ofQ.  with  a  covering  constant  Nc,  i.e.  we  can  color 
{^]iL\  using  at  most  Nc  colors  in  such  way  that  the  subdomains  of  the  same  color 
are  disjoint. 


•  There  exists  a  partition  of  unity  {0,}  satisfying 

^di  =  1,  with  e,  e  C^(fi,),0  <  ^,  <  1  and  \V0,\  <  C/H. 

One  way  of  constructing  {A,}  satisfying  the  above  assumption  is  described  in  Dryja  and 
Widlund  [6].  We  extend  each  fi,  to  a  larger  region  fi,  so  that  CiH,  <  dist{dCli,dQ,)  <  C2Hi. 
We  cut  off  the  part  of  fi^  that  is  outside  Q..  Another  way  of  constructing  { fi,}  is  described 
in  section  4. 

Let  Vo  =  Vq    and  Vk  be  the  bicubic  elements  associated  with  the  triangulations  T^ 

and  T'',  respectively.    Let  V,  =  V^{(li)  =  V^j  n  i/^(fi,).  and  Pv,  :  7/^(0)  -^  V.,  be  the 
orthogonal  projection,  and  let 

TV 

We  have  the  following  additive  Schwarz  algorithm 
Algorithm  3.1.  Find  Uh.  G  V'^  such  that 

(7)  Puh  =  gh, 

with  gh.  =  J2i9i'  where  gi  is  given  by  the  solutions  for  the  following  finite  element  problems 

(8)  a{gi,<f>h)- a{Pv,u,<f>h)  =  f{4>h),     y(l>h  ^  Vi- 

To  find  u/i,  we  first  find  the  right  hand  side  gn  by  solving  (8)  and  we  then  use  the  conjugate 
gradient  method  to  solve  the  system.  In  each  iteration,  we  need  to  compute  Pvh  for  some 
element  v^  G  V^.  This  is  done  in  the  following  steps 

•  Compute  PoVh  by  solving  the  finite  element  problems  K^x  =  6  in  subspace  Vq  . 
Here  K^  is  the  stiffness  matrix  and  dim(Vg^)  =  4A^,  where  N  is  the  number  of 
interior  coarse  grid  points. 

•  Compute  PiV/i  for  each  subdomain  f2,  by  solving  the  finite  element  problems  A'fx  = 
6  in  the  subdomain  A,.  Here  K-^  is  the  stiffness  matrix  and  dim(Vi)  =  4n,,  where 
n,  is  the  number  of  fine  grid  points  inside  Qj. 

Theorem  3.1.   The  iteration  operator  P  of  the  additive  Schwarz  scheme  satisfies 
Cia{u,u)  <  a{Pu,u)  <  C2a{u,u),     Vu  G  Vq. 


It  is  of  course  desirable  to  reduce  the  work  in  each  iteration  without  decreasing  the  rate 
of  convergence.  Let  A^{D)  and  A^{D)  be  the  sets  of  fine  and  coarse  grid  points  in  the  set 
D,  respectively.  Let  </)^  and  $^  be  the  nodal  basis  functions  of  Vq  and  Vq  ,  respectively.  To 
describe  an  alternative  algorithm,  we  also  need  to  introduce  the  following  new  subspaces 

•  The  vertex  space  V^^  =  sp3iTi{4>^^{x)},  for  each  vertex  x^.  For  </>  G  V^^,  T{xk)  = 
supp((^)  is  a  polygon  of  diameter  0{h).  We  note  that  dim(Vj^^)  =  1. 

•  Vi  =  spa.n{<f>^{x),\a\  <  l,Xk  G  A''(Qi)}.  V  is  a  subspace  of  Vi  consisting  of 
functions  with  vanishing  second  derivative  (g^  )  at  all  the  vertices  of  the  elements. 
The  stiffness  matrix  for  Vi  is  a  principal  minor  of  the  stiffness  matrix  for  Vj  and 
dim(V;)  =  f  dim(V;). 


•  Vq    =  span{$f,|a|  <  l,i  G  A'^}.  Vq    is  a  subspace  of  Vq    consisting  of  functions 
with  vanishing  second  derivative  {-q^  )  at  all  the  vertices  of  the  substructures. 
The  stiffness  matrix  for  Vq   is  a  principal  minor  of  the  stiffness  matrix  for  Vq   and 
dini(y^)  =  f  dim(y«). 
Using  these  subspaces,  we  have 

Algorithm  3.2.  Find  Uh  G  V^  by  solving 

(9)  p^'^uH  =  (Py^  +  E  ^v.  +  E  Py^, )"''  =  9h 

with  the  appropriate  right  hand  side  gh- 

Theorem  3.2.    The  iteration  operator  P  for  Algorithm  3.2  satisfies 
Cia(u,  u)  <  aiP^^K,  u)  <  C2a{u,  u),     Vu  G  V^ 


3.2.  Argyris  and  Bell  Elements.  As  in  the  bicubic  element  case,  we  define  two 
triangulations,  the  coarse  triangulation  T^  =  {Jl,}  and  the  fine  triangulation  T''  =  {r,}, 
using  triangular  elements.  The  subregions  fi,  are  defined  similarly.  We  assume  that  all  the 
substructures  and  elements  are  shape  regular  in  the  usual  sense.  We  present  our  algorithms 
for  the  Bell  element.  The  algorithms  for  the  Argyris  element  are  similar. 

Let  V^  and  V^  be  the  space  of  Bell  elements  with  respect  to  T^  and  T^,  respectively. 
V^  and  V^  are  similarly  defined. 

In  general,  the  second  derivatives  of  $  G  Vg  at  the  edge  nodes  x,  have  two  different 
values  except  at  the  vertices  of  the  substructures.  Therefore,  Vg  ^  Vg.  Thus,  a  new  coarse 
space  has  to  be  found.  An  easy  way  of  modifying  Vg  to  achieve  this  goal  is  by  replacing 
the  basis  functions  of  Vg  .  Note  that,  in  a  substructure  ilj,  a  basis  function  $  of  Vg  can 
be  represented  by  the  basis  of  Vg: 

^w=    E    E  ^"(^•)<^n^)+    E     E  ^"(^•)<^"(^)  ^e^^' 

x.€A''(n,)|a|<2  Xi€A''(afij)|a|<2 

which  is  now  replaced  by 

*w=    E    E^-ixi)4>nx)+    E    E*-(^')<(^)- 

x,eA'>(nj)|o|<2  ieA''(afij)  |a|<2 

Note  that  the  second  derivatives  of  $  vanish  at  the  nodes  on  dCti-   We  define  the  coarse 
space  Ug  as 

U§  =  spa.n{9f,\a\<2,ieA"}. 

It  is  easy  to  see  that  we  have  the  inclusion  U^  C  V^.  Let  Vq  =  U§ ,  V;  =  V^  n  ^o(^«)  ^nd 
we  obtain  a  space  decomposition 

N 

f''  =  ek 

t=0 

and 


Algorithm  3.3.  Finduh  G  Vg  such  that 
(10)  Puh  =  {Pvo  +  Pv,  +  ---  +  Pvs  )uh  =  9h 

with  an  appropriate  gh,- 

Theorem  3.3.    The  iteration  operator  P  for  Algorithm  3.3  satisfies 
Cia{u,  u)  <  aiP^^K,  u)  <  C2a{u,  u),     Vu  e  Vj 

Remark  3.1.  In  computing  PoVh,  we  need  to  solve  the  coarse  problem  KyHX  =  6, 
where  K^jh  =  {a('J',,  *j)}  and  '*,,  *j  the  modified  basis  functions.  Thus,  we  need  to 
compute  the  matrix  Kjjh  =  {a{^i,9j)}.  A  standard  way  is  to  use  numerical  integration. 
An  alternative  way  is  to  replace  K^h  by  KyH.    This  is  equivalent  to  using  an  inexact 

solver  Pq  to  replace  Pq.   It  can  be  shown  that  {a($i,$_,)}  and  {a('f, ,$_,)}  are  spectrally 
equivalent.  Thus  we  still  have  an  algorithm  with  uniformly  bounded  condition  number. 

There  are  also  simplified  algorithms  for  the  Argyris  and  Bell  elements  which  are  quite 
similar  to  those  for  the  bicubic  element. 

4.  Multilevel  Methods  for  the  Biharmonic  Problem.  In  this  section,  we  consider 
multilevel  additive  Schwarz  methods  for  the  biharmonic  equation.  Although  all  the  two 
level  algorithms  in  this  paper  can  be  generalized  to  more  than  two  levels,  we  only  consider 
a  special  case,  which  in  matrix  form  corresponds  to  a  multilevel  block  diagonal  scaling 
(MBDS).  For  simplicity,  we  use  the  bicubic  elements. 

We  define  a  sequence  of  nested  rectangular  triangulations  {T'l^j.  We  start  with 
a  coarse  triangulation  T^  =  {t/},_'j,  where  r/  represents  an  individual  rectangle.  The 
successively  finer  triangulations  T'  =  {t/},_!j  are  defined  by  dividing  the  rectangles  of  the 
triangulation  T'~^  into  four  rectangles.  Let  h\  =  diameter(r'), /i/  =  max,/i',  and  h  =  hi,. 
The  level  /  grid  points  are  denoted  by  A',  and  the  basis  functions  by  </''a,i  G  A'.  Here  a 
represents  the  order  of  derivatives. 

Let  V'  =  Vq'  be  the  bicubic  elements  associated  with  T'.  Let  fi'  =  supp{(/)'^}  be  the 
support  of  an  individual  basis  function,  and  let  V-'  =  span{</>'  ^}  be  the  span  of  the  level  / 
basis  functions  at  the  grid  point  i,.  We  note  that  for  the  bicubic  element  dim{V/}  =  4.  On 
each  level,  we  have  an  overlapping  decomposition  of  the  domain 

n  =  u,ii|. 

This  decomposition  satisfies  Assumption  3.1.  We  have  the  space  decomposition 

1=1  i=l 

The  operator  of  the  i-level  additive  Schwarz  algorithm  is  given  by 

L     N,  L     N, 

(=1 i=i     1=1  i=i 

where  Pyi  is  the  a(-,  •)-orthogonal  projection  from  V^  onto  V-'. 
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Algorithm  4.1  (MBDS  Algorithm).  Find  Uh  e  V^  by  solving 
(11)  Puh  =  9h 

with  an  appropriate  right  hand  side  gh- 

In  the  matrix  form,  equation  (11)  can  be  written  as: 

B-^Klx  =  B-H 

where  5"^  =  Dl^  +  ITL-ii^J^liII^.i  +  •••  +  niA'f^II'i.  Here  Ki  is  the  stiffness  matrix 
associated  with  T'  and  Di  -  diag{A'/},  IT;  a  prolongation  operator,  and  Ilf  a  restriction 
operator.  The  MBDS  algorithm  is  a  natural  generalization  of  a  block  diagonal  scaling 
method. 

Theorem  4.1.    The  multilevel  additive  Schwarz  operator  P  satisfies 
C\L~^a{u,u)  <  a{Pu,u)  <  CLa(u,u). 
Thus 

Remark  4.1.  For  second  order  problems,  the  corresponding  algorithm  is  a  multilevel 
diagonal  scaling,  which  is  equivalent  to  the  BPX  algorithm  of  Bramble,  Pasciak  and  Xu  [1]. 
We  note  that  in  the  second  order  case,  we  have  that  the  operator  P  has  a  constant  upper 
bound;  cf.  Zhang  [10].  It  is  also  possible  to  strengthen  the  result  in  Theorem  4.1. 

REFERENCES 

[1]  J.  H.  Bramble,  J.  E.  Pascuk,  and  J.  Xu,  Parallel  multilevel  preconditioners,  Math.  Comp.,  55 
(1990),  pp.  1-22. 

[2]  X.-C.  Cai,  Some  domain  decomposition  algorithms  for  nonselfadjoint  elliptic  and  parabolic  partial 
differential  equations.  Tech.  Rep.  461,  Computer  Science  Department,  Courant  Institute  of  Math- 
ematical Sciences,  September  1989.  Courant  Institute  doctoral  dissertation. 

[3]  ,  An  additive  Schwarz  algorithm  for  nonselfadjoint  elliptic  equations,  in  Third  International  Sym- 
posium on  Domain  Decomposition  Methods  for  Partial  Differential  Equations,  T.  Chan,  R.  Glowin- 
ski,  J.  Periaux,  and  O.  Widlund,  eds.,  SIAM,  PhUadelphia,  PA,  1990. 

[4]  X.-C.  Cai  and  O.  Widlund,  Domain  decomposition  algorithms  for  indefinite  elliptic  problems,  SIAM 
J.  Sci.  Statist.  Comput.,  13  (1992).  To  appear. 

[5]    P.  G.  ClARLET,  The  Finite  Element  Method  for  Elliptic  Problems,  North- Holland,  1978. 

[6]  M.  Dryja  and  O.  B.  Widlund,  An  additive  variant  of  the  Schwarz  alternating  method  for  the  case  of 
many  subregions,  Tech.  Rep.  339,  also  Ultracomputer  Note  131,  Department  of  Computer  Science, 
Courant  Institute,  1987. 

[7]  ,  Some  domain  decomposition  algorithms  for  elliptic  problems,  in  Iterative  Methods  for  Large 

Linear  Systems,  San  Diego,  California,  1989,  Academic  Press,  pp.  273-291.  Proceeding  of  the 
Conference  on  Iterative  Methods  for  Large  Linear  Systems  held  in  Austin,  Texas,  October  19  -  21, 
1988,  to  celebrate  the  sixty-fifth  birthday  of  David  M.  Young,  Jr. 

[8]  ,  Multilevel  additive  methods  for  elliptic  finite  element  problems,  in  Parallel  Algorithms  for  Par- 
tial Differential  Equations,  Proceedings  of  the  Sixth  GAMM-Seminar,  Kiel,  January  19-21,  1990, 
W.  Hackbusch,  ed.,  Braunschweig,  Germany,  1991,  Vieweg  fc  Son. 

[9]   A.  M.  Matsokin  and  S.  V.  NEPOMNYASCraKH,  A  Schwarz  alternating  method  in  a  subspace,  Soviet 
Mathematics,  29(10)  (1985),  pp.  78-84. 
[10]   X.  Zhang,  Multilevel  additive  Schwarz  methods,  tech.  rep.,  Courant  Institute  of  Mathematical  Sciences, 
Department  of  Computer  Science,  1991.  To  appear. 

[11]  ,  Studies  in  Domain  Decomposition:  Multi-level  Methods  and  the  Biharmonic  Dirichlet  Problem, 

PhD  thesis,  Courant  Institute,  New  York  University,  1991.  To  appear. 


NYU  COMPSCI  TR-583 

Zhang,  Xuejun 

Domain  decomposition 
algorithms  for  the 
biharmonic  Dirichlet 


c.l 


_NYU    COMPSCI    TR-583     ^    ]_ 
"zhang,    Xuejun 
_Domain   decomposition 

algorithms   for    the                    ~ 
_  biharmonic   Dirichlet 

DATE  OUE 

BORROWER'S  NAME 

LIBRARY 

N.Y.U.  Courant  Institute  of 

Mathematical  Sciences 

251  Mercer  St. 
New  York,  N.  Y.    10012 


This  book  may  be  kept 

FOURTEEN    DAYS 

A  6ne  will  be  charged  for  each  day  the  book  is  kept  overtiine. 

AYLORD    1*2 

PRINTED   IN   U    S   A 

